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ABSTRACT 


This  study  is  to  understand  the  effects  of  nonlinear 
behavior  of  structures  using  a  finite  element  method. 

The  nonlinear  behavior  equations  are  derived  from 
equations  of  motion  and  constitutive  equations.  The  basic 
theory  is  the  principle  of  virtual  work. 

The  thesis  begins  with  a  comprehensive  formulation  of 
continuum-based  finite  element  theory.  The  theoretical 
portion  concludes  with  the  details  of  both  spatial  and 
temporal  discretization,  including  a  discussion  of 
nonlinearity. 

In  particular,  large  displacement  problems  and  visco¬ 
elastic  problems  remain  a  challenging  engineering  problems 
today.  The  viscoelastic  problems  depend  on  the  relaxation 
function  which  is  the  source  of  material  nonlinearities. 
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A.  ANALYSIS  OF  STRUCTURES 

Analysis  of  structures  subjected  to  large  displacements 
requires  the  inclusion  of  ndnlinear  terms  in  the  strain 
displacement  relations.  In  the  formulation  of  the  problem, 
either  total  Lagrangian  or  Eulerian  approach  is  adopted.  In 
order  to  assure  structural  integrity  for  loads  beyond  linear 
range,  nonlinear  analysis  is  very  much  mandated. 

Another  class  of  problems  is  the  viscoelastic  behavior 
of  structures  -  problems  where  stresses  in  the  components 
vary  with  time  at  a  constant  deformation  state  or  strains 
vary  with  time  at  constant  stress  states. 

These  problems  as  applied  to  simple  one-dimensional 
structures  will  be  addressed  in  this  research  work.  In  what 
follows  is  a  brief  survey  of  literature  in  this  area.  At 
the  end  of  this  chapter,  overview  of  the  thesis  is 
described. 

B.  LITERATURE  REVIEW 

A  method  for  nonlinear  static  and  dynamic  analysis  using 
cable  and  beam  elements  is  given  by  Bergan  et  al  [1]. 
Formulation  allows  for  material  nonlinearities  and  very 
large  deformations.  Gadala  and  Oravas  [2]  present  a 
comparison  of  finite  element  formulations  based  on  the 


Lagrangian,  updated  Lagrangian  and  Eulerian  formulations  for 
nonlinear  mechanics  problems. 

Henriksen  [3]  applies  a  single  integral  constitutive  law 
to  describe  viscoelastic  behavior  in  his  formulation. 
Ghoneim  et  al  [4]  treat  the  total  strain  rate  as  sum  of 
viscoelastic  and  viscoplastic  contributions  in  implementing 
the  constitutive  equations  in  their  formulation,  Pinsky  and 
Kim  [5]  use  variational  approach  for  studying  the  nonlinear 
viscoelastic  shell  behavior  by  representing  the  layered 
shell  as  a  multi-director  field.  The  formulation  includes, 
within  individual  layers,  the  effects  of  transverse  shear 
and  transverse  normal  strain  to  arbitrary  orders  in  the 
layer  thickness  coordinate  [6,  7]. 

Saigel  and  Yang  present  a  quadrilateral  shell  element 
and  describe  implementation  of  elastic-viscoplatic 
constitutive  law  for  large  deformation.  A  large  collection 
of  papers  on  the  application  of  ADINA  may  be  found  in  [8]. 

Response  of  nonlinear  mechanical  systems  to  random 
excitation  is  given  by  To  [9].  Thermomechanical  response  of 
uniaxial  bars  with  thermoviscoplastic  constitutive  laws  is 
discussed  by  Allen  [10].  Large  deformations  of  the 
contacting  bodies  on  the  contact  surfaces  and  solution 
methods  are  given  consideration  in  [11,  12].  Large 
deformations  of  composite  beams  subjected  to  axial  dynamic 
load  is  reported  in  [13]. 

Based  on  Eulerian-Lagrangian  formulation,  frictional 
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contact  analysis  is  discussed  by  Haber  et  al  [14]. 

Halleux  and  Casadei  [6]  present  transient  large  strain 
analysis  of  structures  based  on  biquadratic  elements. 
Polar-decomposition  is  used  in  arriving  at  explicitly 
derived  tangent  stiffness  matrices  for  a  frame  element  by 
Kondoh  and  Atluri  [15].  A  method  based  on  co-rotational 
coordinate  system,  where  the  local  cartesian  coordinate 
system  attached  to  each  element  translate  and  rotate  is 
described  by  Mattiasson  [16]. 

Inelastic  problems  based  on  incremental  virtual  work 
equations  is  given  by  Nagtegaal  [17].  A  constant  strain 
rate  is  assumed  during  an  increment  in  integrating  the 
constitutive  equations. 

Planar  curved  beams  undergoing  large  rotations  is 
analyzed  using  Helinger-Reissner  mixed  variational  principle 
by  Noor  [18].  Steele  discusses  several  nonlinearities 
including  buckling,  creep,  and  nonlinear  elastic  problems  in 
[19].  The  use  of  Lagrangian  formulation  using  an  elasto- 
plastic  model  is  presented  by  Voyiadjis  and  Navaee  [20]. 
Hirashima  et  al  [21]  present  an  approach  for  the  snap- 
through  analysis  of  curved  shell  panels  and  determination  of 
initial  stress  matrix. 

Importance  of  viscoelastic  relaxation  effects  on  a 
pellet  are  investigated  by  Basombrio  [22],  where  it  is 
observed  that  the  duration  of  the  ramp  input  influences  the 
response  considerably.  Viscoplastic  constitutive  models  are 


obtained  explicitly,  without  resorting  to  inversion,  to  be 
directly  implemented  in  a  finite  element  formulation  [23]. 
An  intrinsic  time  scale  that  can  vary  with  time  and  position 
is  used,  together  with  means  to  control  the  numerical 
integration  in  the  approximation  of  constitutive  model  by 
Chambers  and  Becker  [24]. 

Chandra  et  al  [25,  26]  discuss  applications  of  visco¬ 
plastic  models  in  the  metal-forming.  Elasto-plastic  models 
are  used  in  analyzing  unilateral  contact  and  friction 
problems  by  Cheng  et  al  [12].  Chung  et  al  [27]  suggest 
methods  for  varying  parameters  in  guiding  solutions  for  the 
simulation  of  plastic  forming  processes.  Liu  et  al  [28] 
propose  formulation  for  path- dependent  materials  based  on 
arbitrary  Lagrangian  Euler  method  and  include  an  elastic- 
plastic  wave  propagation  problem. 

Rheinboldt  and  Riks  describe  several  algorithms  based  on 
continuation  method  and  mesh  refinements  for  the  solution  of 
nonlinear  finite  element  equations.  Ryu  and  Arora  [29] 
propose  an  approach  to  solve  nonlinear  finite  element 
equations  based  on  substructures  concept. 

Several  papers  pertaining  to  visco-elastic,  visco¬ 
plastic  behavior  and  their  implementation  is  reported  in 
[10,  7,  4,  30,  31]. 
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A.  OVERVIEW 

In  this  chapter  we  consider  a  nonlinear  continuum 
equations  for  a  rod  element.  This  approach  is  based  on  the 
principle  of  virtual  work  [48].  The  resulting  equations  are 
in  convenient  form  for  numerical  solution  via  finite  element 
methods  [49].  This,  of  course,  is  the  conceptual  basis  of 
the  discretized-solid  finite  elements,  introduced  by 
J.T.Oden  et  al  [50]  in  the  elastic  structures  and  nonlinear 
continua. 

B.  NONLINEAR  CONTINUUM  EQUATIONS 

Discuss  a  rod  element  with  the  total  Lagrangian  Method. 
It  may  be  proved  that  internal  virtual  work  is  equal  to 
external  virtual  work.  Namely  a  deformable  system  is  in 
equilibrium  if  and  only  if  the  total  external  virtual  work 
is  equal  to  the  total  internal  virtual  work  for  every 
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virtual  displacement  consistent  with  the  constraints. 

1.  Principle  of  Virtual  Work 

In  this  section  we  prove  that  total  internal  virtual 
work  is  equal  to  total  external  virtual  work.  Consider  the 
concept  of  minimum  potential  energy  [51].  The  condition 
is  satisfied  to  equilibrium  state  and  the  change  of 

potential  with  respect  to  the  displacement 
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consider  the  load- 


must  remain  stationary.  And 
displacement  and  stress-strain  relationship  of  figures 

shown  below. 


Figure  2.1  Load-Displacement  Curve 


Figure  2.2  Stress-Strain  Curve 

In  General,  the  total  potential  energy  is  equal  to  the 

summation  of  strain  energy  and  potential  energy. 

ie,  xp  =  U  +  W  (2.1) 

where,  xp  :  Total  potential  energy 

U  :  Strain  energy  (internal  work) 

W  :  Potential  energy  of  load  system  (external 
work) 
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From  Figure  (2.1), 


w  ,-J 


p  6u 


=  -  2  Pi  Ui 

where,  Pi  :  ith  load 

ui  :  deflection  under  i*h  load 


(2.2) 


From  Figure  (2.2), 

«  =  J  f‘ 


a  6€  d(vol) 


VO  1  o 


...J 


a  €  d(vol) 


VO  1 


(2.3) 


Substituting  the  equation  (2.2)  and  (2.3)  in  equation  (2.1), 
we  get 


,  f.  ■ 


rr  =  J  .5  a  €  d(vol)  -  2  p  u 

P  vol  Fi  i 


(2.4) 


Equation  (2.4)  can  be  represented  by  matrix  form,  the 
following  is  the  matrix  formula. 


-J. . 


*p  =VJ1  -5  {a)  {€)  d(vol)  ~  {  u4  }  {  pt  > 


(2.5) 


where,  (a) 
(€} 
{ui  } 
(Pi  ) 


a  vector  of  stresses 
a  corresponding  vector  of  strains 
a  vector  of  elemental  nodal  displacements 
a  vector  of  elemental  nodel  forces 


Let  us  dafine  [B] ,  a  matrix  relating  strains  and 
displacements, 


ie,  {€)  =  [B]  {u}  (2.6) 

and  define  [D],  a  matrix  relating  stresses  and  strains, 
ie,  {a}  =  [D]  {€>  (2.7) 


r 


Substituting  equation  (2.6)  &  (2.7)  in  equation  (2.5',  we 

get 


f  T  T  T 

*P  =  J  {Ui  }  {B}  [D]  [B]  d(vol)  {Ui  }  -  {Ui  }  (Pi  >  (2.8) 

VO  1 

Now  we  consider  the  minimum  potential  energy,  the  condition 
is  as  follows: 

dJtp 

-  =  0  (2.9) 

d{ux  } 

Taking  derivatives  with  respect  to  {ui  }  in  equation  (2.8) 
and  applying  to  equation  (2.9),  we  get 

(Pi  >  =  J  [B]T  CD]  [B]  d(vol)  {ui  >  (2.10) 

VO  1 

From  relating  loads  and  displacements,  we  can  get  the 
elemental  stiffness  matrix  [k] . 

ie,  [K]  =  J  [B]  [D]  [B]  d(vol)  (2.11) 

vol 

Summing  all  the  elements  of  the  entire  structure,  equation 
(2.10)  becomes 

{P>  =  [K]  {u>  (2.12) 

where,  {P>  =  2  {Pi  > 

[K]  =  2  [K] 

[u]  -2  {ui  } 

2.  One  Dimensional  Finite  Element 

Consider  a  shape  function  which  is  linearly  function 
with  geometric  boundary  condition  with  one  dimensional 
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finite  element.  We  define  that  the  shape  function,  *«.  ii 


=  aM  +  bM  X 


In  order  to  derive  the  constants  in  the  shape  function  for  a 
rod  element  with  linear  displacement  field,  use  the 
geometric  boundary  condition. 

We  may  obtain  the  constants  as  follows: 

X2 


Lo 


1 


Lo 

xi 

Lo 

1 

Lo 


(2.13) 


The  relationship  between  displacement  and  shape  function  is 
given  below 

M 


u  =  u*  (H=l;  M=2 ,  .  .  .  .  ) 

=  'l'l(x)  Uj  (t)  +  ^(x)  u2  (t) 


(2.14) 


Another  consideration  is  the  viscoelastic  materials.  For 
viscoelastic  materials,  the  stresses  may  be  expressed  in 
terms  of  relaxation  functions  as  given  by 


S( 


t)=  I  0(t  -  £  )  S*  ({)  d£ 
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=  s  (t) 


I  ( t) 


(2.15) 


e  e 

where,  S  (t)=S  [r(t)] 


I( 


:  elastic  response 

r  g  ( 1  -  £  >  t  t 

=  J - TT- .  Se  (|)  df 

°  <9? 


G(t):  relaxation  function. 

Assume  the  element  is  in  equilibrium  and  quasistatic  state. 
The  motion  equation  for  a  one  dimensional  finite  element  is 
given  by  J . T . ODEN  [ref. 6], 


J 


Vo 


(  1  ^  u  )  dV  =  P  (u  ,t) 

N ,  X  M,  K  M  o  NM 


where,  Vo  =  AoLo  :  undeformed  volume 

P  :  total  generalized  force  at  node  N 

M 

On  substituting  for  above  shape  functions  and  their 

derivatives,  we  may  obtain 

( *-*M  ’  ^  ^  C  S*  “  I  ]  (1  *  U^  )  “  P||  ( >  t ) 

=  0  (2.16) 

C.  CONSISTENT  NODAL  FORCES 

Consider  a  rod  element  in  two  configurations:  undeformed 
and  deformed. 

Let  the  distributed  forces  per  unit  length  in  these 
configurations  be  q  and  q*  respectively.  The  consistent 
load  vector  depends  on  both  the  magnitude  of  the  load  and 
the  deformed  state. 
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One  dimensional  rod  element  under  distributed  load  is 


showed  as  following  figures. 


%*C*,  ■£) 


■ 

Va 

m 

2  pj 

u 


where,  Lo :  undeformed  rod  element  length 
L  :  deformed  rod  element  length 

q(x,t)  is  distributed  load  per  unit  undeformed  length 
q*(y,t)  is  distributed  load  per  unit  deformed  length 

Figure  2.3  One -Dimensional  Rod  Element  under  Distributed 
Load . 


Consider  the  relationship  between  the  deformed  and 
undef  ortned  coordinate.  It  may  be  written  as 

y(x,t)  -  x  ♦  u(x,t)  (2.17) 


where  x:  undeformed  coordinate  of  a  point 
y:  deformed  coordinate  of  a  point 
u:  displacement  of  a  point  within  element 


The  consistent  nodal  forces  vector  is  given  by 

P  =  /vpf-NdV  -j [/AfcJdV. 

In  index  notation,  nodal  forces  may^e  written  as 


(2.18) 


(2.19) 


with  N  =  1,  2 

body  forces  per  unit  mass  deformed 
body  forces  per  unit  mass  undeformed 
deformed  density 
undeformed  density 

Relating  q*(y,t)  and  q(x,t)  is  given  by 
q( x , t )  =  \  q*(y, t) 
where,  \  =  dy/dx  :  stretch 

Consider  the  undeformed  body  force  per  unit  mass, 
may  write 

ft-  _  y*da.  _  %  _  x  V 

r* '  dm  "  p.  A.  '  p.  A. 

Using  equation  (2.17) 

y  =  x  ♦  u 

M  M 


where  F 
Fo 

P 


dy  T 

-  =  1  +  U/’  U 

dx  M,  I  M 


The  undeformed  body  force  per  unit  mass  is  then 

q* 

0  o  Ao 


ft- 


(  1  H,  ■  UM  ) 


The  consistent  nodal  forces  take  the  following  form, 

PN  =  ^  H^Sdx 


Substituting  the  shape  function  in  equation  (2.25), 
get 


(2.20) 

then  we 

(2.21) 

(2.22) 

(2.23) 

(2.24) 

(2.25) 

we  may 
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(2.26) 


PLo 

P„  =  (1  +  b^^)  J  q*(aN  +  bH  x)  dx 

Now  we  consider  some  special  cases  of  equation  (2.26). 
One  case  is  the  uniformly  distributed  loads. 


Figure  2.4  Rod  Element  under  Uniformly  Distributed  Loads 


The  uniformly  distributed  loads  are  given 
q*(y,t)  =  qo  h(t) 

Substituting  the  loads  in  equation  (2.26),  we  may  get  the 
consistent  nodal  forces. 


(  PM  >  = 


where,  h(t)  = 


r  P.1  qo  h(t)  f  1* 

■[  J-— "•♦'--’14 


(2.27) 


du 


dx 


Another  case  is  the  linearly  distributed  loads. 

k 


Figure  2.6  Rod  Element  under  Linearly  Distributed  Loads 
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Linearly  distributed  loads  are  given 


q(x,t)  =  q  (t)  +  [q  (t)  -  q  ( 1 )  ^  j 

Substituting  the  loads  in  equation  (2.26)  we  may  get  the 
consistent  nodal  forces. 


D.  SOLUTION  METHODS 

In  this  section  we  derive  the  equilibrium  equations  of 
incremental  method,  Full  Newton-Raphson,  and  Modified 
Newton-Raphson  Method. 

It  is  based  on  total  Lagrangian  Method  and  use  the 
Taylor  Serier  at  a  given  equilibrium  point. 

The  element  stiffness  matrices  and  consistent  load 
vectors  are  assembled  to  form  global  matrices  in  the  usual 
manner.  The  prescribed  displacements  are  implemented  by 
replacing  the  diagonal  element  of  the  corresponding  row  by 
unity  and  the  rest  of  the  element  of  the  rows  and  columns  by 
zeros  [49]. 

1 .  Incremental  Method  (INC) 

Consider  the  incremental  process,  we  can  rewrite  the 
equation  of  load-displacement  relationship, 

[K]»  =  (A  P)i  ♦  { A  Pc  }t  (2.29) 

where,  i  :  the  current  configuration  coordinates. 

(  A  Pc  }i  :  the  force  unbalance 
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{AP}i*i  :  an  increment  in  the  external  load 
The  method  is  described  by  the  Figure  2.6 


where,  ABC  :  The  algorithm  of  equation  (2.29) 

AB’  C’  :  The  algorithm  of  purely  incremental 
process 

Figure  2.6.  The  Algorithm  of  Incremental  process 

The  approximate  values  of  A,  E,  C,  is  computed  fror.,  p  =f(u). 
To  improve  these  values,  take  the  iteration  process  at  each 
load  level  with  {AP}»*i=  0,  as  is  Figure  2.7. 

The  iteration  processes  are  two  ways,  one  is  few 
load  levels  with  many  corrective  iterations  in  each,  or 
another  is  many  ioad  ievels  with  few  iterations  ir.  each. 
This  xs  sometimes  called  incremental  with  one-step  Newton- 
Fcaphson  method" . 

2.  Full  Newton -Raphson  Method  ( FNR ) 

Consider  the  load  versus  displacement  function  P  = 


f(u).  Take  the  Taylor  series  of  P  =  f(u) 


f(UA+  A,)  -  +  (^-)A  A, 


(2.30) 


dP 

where,  -  =  K  =  tangent  stiffness  at  A 

du  A 

A 

f  (  U  )  =  P 

A  A 

f (u  +  A  )  =  P 

A  l  B 

We  look  for  A.such  that  f(u  ♦  A  )  =  P  ,  so 

1  A  1  B 

K  (A)  =  P  "  P  (2.31) 

A  I  B  B 

The  consequent  steps  are  explained  by  the  Figure  2.7. 


Figure  2.7.  The  Newton  Kaphson  solution  of  the  equation 
P  =  f (u) 

The  above  algorithms  are  below; 

.  update  displacements  :  ui  =  Ua  ♦  A , 

.  use  ui  to  get  stiffness  ki  and  load  Pi 
.  find  next  increment  from  ki  (A*)  =  Pb  -  Pi 

Eventually  ua  A,  *■  Aa  ♦  .  =  ub  to  a  close 

approximation . 

Then,  the  load  can  be  increased  from  Pb  to  Pc,  and  the 
process  can  be  started  again.  This  solution  algorithm  is 


known  as  the  Newton-Raphson  Method. 

3.  Modified  Newton-Raphson  Method  (MNR) 

The  Newton-Raphson  method  requires  that  tangent 
stiffness  matrix  [k]  be  generated  and  reduced  for  equation 
solving  in  every  iteration  cycle.  This  process  is  not 
economic.  An  alternative  is  modified  Newton-Raphson 

iteration,  also  called  "constant  stiffness  iteration".  Here 
the  same  tangent  stiffness  matrix  is  used  for  several 
iterations  cycles.  The  MNR  is  described  by  the  Figure  2.8 

& 


Figure  2.8  The  MNR  Iteration 

As  compared  with  FNR,  the  MNR  requires  more  iterations, 
but  each  is  done  more  quickly  because  [K]  is  formed  and 
reduced  in  only  the  first  iteration. 

E.  ANALYTIC  SOLUTIONS 

At  first  we  define  the  concept  of  stresses  [52]  Eulerian 


stress,  o  ,  is  defined  as  follows: 


N(y,t) 

a  (y,t)  =  - 

A(y) 


(2.32) 


and  is  called  True  stress.  Lagrangian  stress,  T,  is  defined 
as  follows: 


N(y , t ) 

T(x.t)  =  -  (2.33) 

Ao  (x) 

and  is  called  Engineering  stress. 

Kirchloff  stress,  S,  is  defined  as  follows: 

N(y, t) 

S(x.t)  =  — -  (2.34) 

\(  x ,  t )  Ao  ( x ) 

The  relationship  between  strain  r  ,  and  stretch  is  as 

follows : 


r  =  S  (  -  1)  (2.35) 

du 

A  =  1  +  -  (2.36) 

dx 

where,  N  :  Applied  force 
A  :  deformed  area 


For  the  special  case  of  small  strain,  both  u,  x  and  change 
in  cross-sectional  area  are  small.  In  this  case,  we  can 
write 

r  «*  u  _  =  € 

I  * 

o  "•  T  ••  S 

where,  G  is  engineering  strain. 

Thus  in  small  deformation  case,  the  Lagrangian  stress  T 
and  engineering  strain  G  can  be  used  to  replace  the 
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and  engineering  strain  €  can  be  used  to  replace  the 
Kirchhoff  stress  S  and  Lagrangian  strain  r  . 

1.  Elastic  Materials 

Consider  a  one-dimensional  rod  with  varying 
cross-  sectional  area,  Ao  =  Ao(x),  in  equilibrium  under 
distributed  external  loading. 


(a)  (b) 

where,  f  :  external  distributed  load  per  unit  length 

T  :  Lagrangian  stress  within  rod 
Figure  2.9.  One-Dimensional  Rod  under  Axial  Loading 


For  static  problems  and  elastic  materials,  both  f  and  T  are 
independent  of  time, 

■i  r  *  -  £  (  \  J  **  - 

*  ^  *  *  ~  *  \  A  $  A  ~  «  \  •»  f 

From  figure  2.9(b),  the  static  equilibrium  of  forces  in  the 
x-direction  i3 

£  F(x) » -TAo  ♦  f  Ax  ♦  (T  +  (dT/dx)  ax  )  (Ao  ♦  (dAo  /  dx)  ax  )  *  0  ( 3  7  J 

Expanding  in  Taylor  series, 


i 

i 
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£  iT  i£t  =  0 

i  +'7%  +p*4x  T<Tx  4x 


(2.38) 


For  a  continuous  rod, 


i  +  ^CTA,')  = 


in  the  limit, 


(2.39) 


Substitution  of  equations  (2.34),  (2.35),  and  (2.36)  in  the 
above  (2.40) 


where,  u(x)  :  displacement  w.r.t.  x 

u, *  :  first  derivative  of  u(x) 

u, x i  :  second  derivative  of  u(x) 

Note  that  equation  (2.40)  gives  the  values  of  distributed 
force  f  (per  unit  length)  required  to  produce  the 
displacement  field  u(x)  of  a  one-dimensional  rod  once  the 
boundary  conditions,  the  stress-strain  law,  and  varying  area 
function  Ao ( x )  have  been  specified.  Consider  the  rod  with 
constant  area  Ao ,  equation  (2.40)  becomes 


f(x)  a  -Ao  u,xx  1(1  ♦  u,x )  ( dS*  /dr )  ♦  S#] 


(2.41) 


Some  special  cases  of  equation  (2.41)  is  considered.  One 
case  is  the  specified  displacements. 
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where,  u(x)  =  ao  +  a,  x  +  a*  x2 
u,  x  =  a,  +  2  at  x 
u,  ii  -  2  a. 

Figure  2.10  Test  Model  of  "Inverse  Solution' 


Apply  the  boundary  contition  u,x(x)  =0  at  x  =  Lo 


u  ( x )  =  a  ,  x  - - x2 

2Lo 


(2.42) 


If  a,  is  specified,  we  may  get  the 


exact  solution  u(x)  and 


the  corresponding  required  distributed  load  f(x). 
a.  Linear  displacements  field: 
u(x)  =  a, x 
u,  x  =  a, 
u,  xx  =  0 


|a 

$ 

'.**<« 


From  equation  (2.42),  we  may  get 
f (x)  =  0 

Note  that  a  linear  displacement  field  cannot  be  produced 
with  distributed  load. 


b.  Quadratic  displacement  field: 


u(x) 


i,  x  -  (  in  /2Lo  )  x2 


u,  *  =  a,  -  y  (  a  i  /  Lo  )x 

u ,  x  x  =  -  ( a  |  /  Lo  ) 

-  Consider  the  linear  materials,  S*  =  E  r 

We  may  get  the  distributed  forces  from  equation  (2.41) 
f  ( x )  s  -  EAuii,  m  [(1  +  u,  x  )  +  ( u,  *  ♦  ^  u2  ,  *  )  ]  (2.43) 

-  -  Consider  the  nonlinear  materials, 

Se  =  El  r  +  E2  r2  (quadratic  form) 

From  equation  (2.41),  we  may  get 

f(x)  =  -  Aou.xx  {  (1  ♦  u.  x  )  [El  ♦  2  E2  ( u,  *  ♦  \  u2 ,  *  )  ] 

+  El  (u.«  f  ^U2,*)  +  E2  (u.  *  •*-  «au2.,)]  >  (2.44) 

Another  case  is  the  prescribed  loads. 


Figure  2.11  Test  Model  of  Rod  with  Prescribed  End  Element. 


From  equation  (2.39),  since  no  distributed  loads  exist, 


(TAo)  =  0  with  constant  area. 


We  may  get 


S  =  constant 


(2.45) 


»  .*  ,1  «  i  | 


(2.46) 


The  applied  load  is  given 
P  —  Ao  \  ii 
-  Linear  material; 

r  £r  -  E(  u,  »  ♦  >2  u- .  *  ) 


From  equation  (2.46  i,  we  may  get 

P  =  EAo  (u,«  ♦  u2,«  M  uJ,«)  (2.47) 

-  Nonlinear  material; 

=  Fir  *-  E2  r2  (quadratic  form) 

From  equation  (2.46)  we  may  get 

P  -  Ao  (1  ♦  u,  >  )  f  El  (ut  i  ♦  >5  u:  ,i  )  +  E2  (u,  s  ♦  u2  , >  )a  ] 

(2.48) 


2.  igyasilinsjar  Vi^ooela 

The  Viscoelastic  materials  has  the  more  complicated 
problems  than  the  elastic  materials,  since  the  external 
loads  and  displacements  are  time-dependent. 


U.t*> 

*•(*>  - 

< -  - 

Y 

K2 - 


whe re , 


Figure 


Ao ,  Lo  :  Underformed  Area  and  Length 
uo(t)  :  time-varying  left-end  displacement 
uLo(t)  :  time-varying  right-end  displacement 
Fo  (t)  :  time-varying  left-end  force 

FLo(t)  :  time-varying  right-end  force 

2.12  One-Dimensional  Rod  Subjected  to  Time-varying 
End  Loads 


For  quasilinear  viscoelastic  materials, 


S(t)=S*U)  -  f 

Jo  <* 


(2.15) 


The  Lagrangian  stress  is  constant  throughout  the  rod  since 
the  distributed  loads  is  not  applied.  From  equation  (2.33), 
the  stress  is  time  variable  only. 

Ti x, t )  -  g( t )  =  \( t )  d( t )  (2.49) 

And  the  internal  forces  prudueed  with  the  rod  is 

F<  t)  -  TAo  =  Aog(t )  (2.50) 

Substituting  equation  (2  15)  m  equation  (2.49)  yields 
g(t)  -  \ ( t )  Sit) 

=  (1  +  «,f)[s*W-  |0~Sp^S#(r)<*|  (2.51) 

The  right  hand  side  of  equation  (2.51)  is  a  time  function 
on  1  y . 

Therefore  g(t)  must  also  be  a  time  function, 
du 

- —  =  u,« (x.t)  =  h(t)  (2.52) 

dx 

Integrating  w.r.t.  x  at  equation  (2.52) 

u  =  u ( x , t )  =  xh(t)  ♦  C(t)  (2.53) 

where,  C  :  integrating  function  to  be  determined  by  the 
boundary  condition  at  one  end. 

Assume  that 

U  (  O  ,  t  )  -  Uo  (  t  )  =  C  (  t  )  , 
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t  lie  re  fore 


u{  x,  t )  =  xn(  t  >  *■  uo  ( t ) . 
u,«  (x,t)  -  hit) 


(2.54) 


substituting  equation  (2.54)  m  equation  (2.51)  yields 


?  i  *  k(*)K  s*{  kill  ♦  i  U*  {k<s)4ft$5l 

1  5  (2.55) 

=  3<c» 

on the  G(t)  and  S* { r )  have  been  specified,  the  equation 
i,:  55)  can  be  used  with  the  boundary  conditions  to  find  the 
exact  solution  in  some  cases.  Consider  the  3  element  solid 
with  linear  material  property.  The  prescribed  displacements 
boundary  conditions  at  both  ends  is  described  in  Figure 
2.13. 


>  X 


BC.  '3  :  u  ( o ,  t )  =  uo(t)  :  0 

u  (Lo  ,t)=  uLo(t)  (2.56) 

figure  2.13  Test  Model  of  a  Rod  with  Prescribed  Time  Varying 
Displacement  at  the  End. 


The  prescribed  displacements  are  linear, 
u(x,t)  -  (x/Lo  )  u,  (t) 

**  O 

u.  «  (X,t)  =  Ut  (  t  )  /  In  =  h(t) 


(2.57) 

(2.58) 


where,  u,  (t)  =  Loh(t) 

“  O 

Assume  linear  material  property,  S*  =  Er 
And  for  a  3  element  rod,  Git) is  defined  below 

G(t)  =  ( 1  -  a  )  ♦  e-Ot  i2.59) 

where,  a  ,  0  ,  constants  for  material  properties. 

Note  that  G(o)  -  1  at  t  =  0 ,  G ( t )  =  1  -  a  at  0  =  0 

No  relaxation  occurs  and  the  material  exhibits  elastic 
response . 

Specifying  the  prescribed  displacement  function  u  .  (t). 

L  o 

a  Ramp  Function 

At  x  =  Lo,  B.  C.’S  :  u(o,t)  =  0 

u(Lo,t)=.u  (t)  =  at  (a:  constant) 

u  O 

From  equation  (2.57)  and  2.58),  we  may  get 

u(  x,  t )  =  ( a/Lo  )  x  t  (2.60) 

u,  *  (  x  ,  t )  =  (  a/Lo  )  t 

Substituting  equations  (2.57)  -  (2.59)  in  equation  (2.55) 

We  may  get  the  reaction  force,  F. 

EAo  at  a^t*  a  a  a  -0t 

F  ;  -  [  1  ♦  - ]  [{at*  - ) - {l-( - )}e 

Lo  Lo  2Lo  (3  0Lo 

a  a  a 

- {(0t  -  l)  -  (02t*  -  2(3t  ♦  2))]  (2.61) 

0  2  0Lo 

Equation  (2.61)  describes  the  reaction  force  produced  at  the 
fixed  end  for  the  prescribed  displacement 

u  r  ( t )  =  at 

L  O 
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b.  Step  Function 


(Relaxation  Test) 


At  x  =  Lo,  B.C.’S  ;  u(o,t)  =  0 
u(  Lo  ,  t )  =  ut^l  ( t ) 

where,  l(t)  •  unit  step  function 
From  equation  (2.57)  and  (2.58),  we  may  get 

ULo 

Ui X, t )  =  -  x  1 ( t )  (2.62) 

Lo 

UL  o 

U.  «  -  -  1  ( t)  (2.63) 

Lo 


Substituting  equations  (2.57)  -  (2.59)  and  (2.62)  in 
equation  (2.55),  we  may  get  the  reaction  force,  F. 


EAo  ULo  UL  o 

F(t)  =  ( - )  ULo  (1  ♦  -  )  G(t)  [1  ( - 


Lo  2Lo 

Consider  the  Kirchhoff  stresses. 
From  equation  (2.34) 

ULo  ULo 

S  =  E  ( - )  (  1  ♦  -  )  G(t) 


)  1  ( t )  ] 

(2.64) 


Lo 


2Lo 


E  rfc  G(t) 


(2.65) 


where,  r,  - 

a  prescribed 
end . 


ULo  ULo 

-  (  i  ♦  - )  :  Lagrangian  strain  due  to 

Lo  2Lo 


step  displacement  function  u.  l(t)  at  the  rod 
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III.  PROGRAM  IMPLEMENTATION 


A.  OVERVIEW 

We  are  interested  in  static  analysis.  For  static 
analysis,  load-displacement  relationship  is  as  following: 
iP}  =  [k]  {u} 

where  { P }  :  external  loads  vector 

[k]  :  constrained  stiffness  matrix 
{u}  :  nodal  displacements  vector 
If  we  know  {P}  or  {u},  then  we  can  analyze  the  static 
problems  by  using  elementary  elastic  theory.  In  general, 
the  static  analysis  steps  are  shown  as  follows: 

Step  1  :  Enter  geometry  of  material  properties,  loads, 
et  al . 

Calculate  element  stiffness  matrices. 

Determine  appropriate  part  of  [k]  and  {P>. 

If  [k]  or  {P>  is  appropriate,  then  repeat  step  2 
through  step  3 . 

If  [k]  or  {P}  is  not  appropriate,  then  apply 
constraints . 

Step  6  :  Determine  {u}. 

Step  7  :  Regenerate  elements. 

Step  8  :  Calculate  element  stresses. 


Step  2 
Step  3 
Step  4 

Step  5 
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B.  SOLUTION  PROCEDURES 

In  nonlinear  structural  problems,  to  find  the  solutions 
we  use  the  INC  and  FNR  and  MNR.  In  fact,  it  is  hard  to  know 
in  advance  what  increments  of  loads  or  iteration  should  be 
used  to  obtain  a  good  approximation  to  the  solution.  The 
interactive  method  was  used  in  determining  the  strategy  for 
solution  in  program  application.  The  strategy  subroutines 
were  developed  for  elastic  or  visco-elastic  time-dependent 
structural  problems.  Here  the  initial  state  for  incremental 
loading  process  was  selected  as  a  null  displacement  vector. 
Next  we  consider  finite  element  method  for  elastic  (FEMEL) 
and  viscoelastic  (FEMEV)  programs .  IN  both  programs, 

equilibrium  was  checked  at  each  stage  for  solutions  obtained 
from  other  INC  or  FNR  and  MNR. 

1  Algorithm  Implementation 

Consider  the  algorithm  of  Newton-Raphson  Method 
based  on  total  Lagrangian  Method.  Full  Newton-Raphson 

Algorithm  is  as  follows: 

Step  1.  Establish  local  coordinates  X-Y  for  the  element  at 
hand. 

Step  2.  Generate  the  element  stiffness  matrix  [k]  in  local 
coordinates  so  that  it  operates  on  the  local  D.O.F. 
Step  3.  Transform  [k]  to  global  coordinates  so  that  it 
operates  on  global  D.O.F. 

Step  4.  Repeat  step  1  through  3  until  all  elements  have  been 
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Step  5.  Compute  displacements. 

Step  6.  Compute  internal  forces,  [k]{u}. 

Step  7.  Solve  [k]  {Au}  =  (A?} 
where,  uN+1  =  u„+  A  uk 

Step  8.  Test  for  convergence. 

If  not  satisfied,  return  to  step  1. 

Where,  convergence  =  ((ARtAR)  /  (AptAp))**0.5  < 
0.001 


,iS  r,’r 


,v  Yr  ?r 


am 


I  v  .  NUUK  MCA  I.  EXAMPLES 


We  consider  two  numerical  examples  for  elastic  and 
viscoelastic  materials. 

A.  ELASTIC  MATERIAL 

1.  Prescribed  concentrated _ forces . 

Specify  the  rod  element. 


* 


a 


tc  -■  22 . 5000  lbs 
A  =0.25  [ ini ] 

L  =10  [in] 

E  =  lOe  [lbs/ins] 

Consider  the  stiffness  matrix  of  rod  element,  [k] .  For 
the  rod  element, 

EA  r  Ki  i  ki  2  -l 

[K]  = -  (4.1) 

L  Lk2i  k22  J 


Consider  the  boundary  conditions, 


P2  = 

22,5000 

at  x 

=  L 

From 

loads  - 

displc 

icemcnts  relationship, 

EA  r 

■  1  -i. 

Ur 

.  jp‘l 

L  1 

-i  J 

,U2. 

LP2. 

Solve 

for  the 

displacements,  U2 , 

Hence 

,  we  may 

get 

LP2 

U2  = 

— 

9.0  ( in ) 

EA 

Prescribed  di splacements 


X 


Consider  the  same  rod  element  before,  and  specified 
the  quadratic  displacements,  u(x), 
u(x)  =  0.02  x  -  0.001  x  2 

The  first  derivative  is 
du 

-  =  0.02  -  0.002  x 

dx 

The  second  derivative  is 
d-u 

-  =  0.002 

dx2 
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From  equation  (  2.41  )  with  the  constant  area,  we  may 

getthe  distributed  loads  for  any  x. 


d2  u  du  dse 

f(x)  =  -  Ao  -  [(1  +  - ) 2  - 


dx2 

3.0  *  10-3  x2 


-  +  S«] 

dx  dr 

-  3.06  x  +  530 


(4.2) 
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COORDINATES 


B.  VISCOELASTIC  MATERIAL 

Consider  the  rod  element  below, 


F.cfc) 


U,C*> 


X 


Eoundary  conditions;  u(0,t)  =  uo  (t)  =0 

u  (  Lo  ,  t )  -  u  (t)  -  u  l(t) 
t  o  L  0 

Lagrangian  Stress;  T{x,t)  =  g(t)  =  Aft)  S(t) 
Internal  Force;  F(t)  =  T  Ao  g(t) 

P 


Kirchhoff  stress;  S  = 


(4.3) 


X  (  X  ,  t )  Ac  ( X  ) 

For  3  element  solids,  the  relaxation  function  is 
lift)  :  (1  -  a)  ♦  a  e-pt  (4.4) 

where,  Oi  and  13  are  constants  of  the  material 
properties . 

Assume  unit  step  function  for  prescribed  displacements. 

At  that  time,  the  boundary  condition  is 
ufo,  t )  =  0 

u ( Lo  ,  t )  =  u,  ( t )  =  u,  lft) 

“  O  u  o 

where,  lft):  unit  stop  function 
We  may  get  the  displacements  function  with  unit  step 
function . 


ufx ,  t )  =  -  ULo  1 ( t ) 

Lo 


(4.5) 


4  5 
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And  the  Kirchhoff  stress  is  specified 


S  =  E  r 


( 


Substituting  equations  (4.4)  -  (4.7)  in  equation  ( 

EAo  ul  o  UL  o 

F(t)  =  -  Ulo  1  +  -  G(t)  1  +  -  l(t)  (< 

Lo  2  Lo  Lo 

Hence,  we  may  get  the  Kirchhoff  stress,  S, 


S  =  rr 


G(t) 


where,  rLo 


UL  o  UL  o 

-  1  +  -  :  Lagrangian  strain 


Lo 


2Lo 


due  to  the  prescribed  step  displacement  function 
u  1  ( t )  at  the  rod  end . 

Equation  (4.9)  describes  the  relaxation  cases, 
equation  (4.8),  we  may  get  the  internal  forces. 


.7) 

4.3), 

-8) 

•  9) 

From 


1 - w - V - 1 - I - 1 - 1 - 1 - 

ri  n  n  n  iq  ro  ro  ro  i 

NOLDNfU  N0LLVXVT3H 


This  study  is  directed  towards  understanding  the 
nonlinear  behavior  of  structures. 

The  formulation  is  based  on  the  principle  of  virtual 
work  We  consider  the  one  dimensional  nonlinear  continuum- 
based  rod  element  that  includes  large  displacements  and 
viscoelastic  effects.  To  get  the  solution  we  use  the 
incremental  process  with  full  Newton-Raphson  and  modified 
Newton-Raphson  Method.  Numerical  solutions  agree  well  with 
exact  solutions. 

Further  studies  may  be  extended  to  two  dimensional 
structures,  such  as  laminated  composite  plates 
Consideration  of  the  viscoelastic  effects  of  composite 
panels  at  high  temperatures  is  another  possible  area  of 


study . 


An  experimental  study  to  verify  the  results  of  the 


present  study  is  also  recommended. 
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